Introducción al análisis de secuencias NGS - Alineamiento
Curso Genética y Genómica en Producción Animal
Dr. José A. Gallardo, Dr(c) Margarita Rivera.
Pontificia Universidad Católica de Valparaísojose.gallardo@pucv.cl
17 June 2021
INTRODUCCIÓN
El paso de alineación o mapeo de secuencias cortas reads producto de la secuenciación NGS a un genoma de referencia se considera una parte integral de los análisis de genomas (Pham-Quoc, Kieu-Do, & Thinh, 2019).
El problema del mapeo de reads consiste en ubicar o alinear dentro de un genoma de referencia previamente secuenciado los millones de reads para luego ensamblarlas a un nuevo genoma (Bak, Bodziony, Migdałek, Pareek, & Żukowski§, 2020) o en el caso de que no exista un genoma de referencia realizar un ensamble de novo.
Esta carga de trabajo se la realiza con diferentes enfoques de alineación, sin embargo, al igual que cualquier software, cada enfoque presenta diferentes formas de ejecución y consecuentemente de análisis.
Para realizar esta tarea se han desarrollado varios programas como lo son, BWA (Li & Durbin, 2009), Bowtie 2 (JHU, 2016) y otros, cuya función principal es hacer coincidir secuencias cortas provenientes de NGS con secuencias de referencia (Bak, Bodziony, Migdałek, Pareek, & Żukowski§, 2020). Actualmente, BWA es un software muy popular para fines de alineamiento, debido a su compatibilidad con datos de Illumina y la integración de distintos algoritmos de alineamiento con mayor rendimiento para reads de distintas longitudes (Li & Durbin, 2009).
El resultado del mapeo o alineamiento se guarda en un archivo d etexto plano con formato SAM (Sequence Alignment Map) o su equivalente comprimido en binario (BAM). Este archivo contiene de forma codificada toda la información de dónde y cómo cada read se alineo al genoma referencia. De el se pueden luego extraer información de variantes para análisis poblacional o estudios de asociación genómica o información de conteo de reads en regiones codificantes para en análisis de expresión de genes.
FORMATO SAM
El formato SAM consta de un encabezado que comienza con el símbolo @ y una sección de alineamiento que contiene la información de cada uno de los reads que alineo al genoma de referencia. Para conocer más acerca del formato SAM puedes revisar el paper que describe el formato o simplmente en wikipedia. Los archivos SAM se pueden analizar y editar con el software SAMtools. A continuación se muestra un ejemplo del formato SAM con los 11 campos obligatorios del alineamiento y 1 campo opcional.
Formato SAM
SOFTWARES DE BIOINFORMÁTICA
BWA
Es un algoritmo de alineación para alinear lecturas de secuencia o secuencias de consulta largas con un genoma de referencia grande como el humano. Elige automáticamente entre alineaciones locales y de extremo a extremo y admite lecturas de extremos emparejados (Li, 2013).
Puedes visualizar las opciones y más información de bwa Aquí
Samtools
Es un conjunto de utilidades que manipulan alineaciones en los formatos SAM (Sequence Alignment Map), BAM y CRAM. Convierte entre los formatos, clasifica, fusiona e indexa, y puede recuperar lecturas en cualquier región rápidamente, puedes tener más información sobre este programa Aquí
Objetivos del trabajo práctico
- Realizar el alineamiento de una muestra en formato .fastq a un genoma de referencia.
Origen de las muestras
Para facilitar el proceso solo se utilizará una muestra la que provieme de la base de datos SRA del NCBI y corresponden a lecturas crudas del salmón del Atlántico Salmo salar en formato fastq, obtenidas por secuenciación de extremos emparejados con un secuenciador Illumina HiSeq2000.
La muestra se alineará al genoma de la mitocondria del salmón del Atlántico puesto que dado que su pequeño tamaño, el alineamiento no debería tardar más de unos pocos minutos.
Etapas del alineamiento.
Descargar Genoma de referencia de la mitocondria e indexar.
Alinear muestra con genoma de referencia.
Analizar el resultado del alineamiento.
Actividades
1 - Conectar a servidor Pomeo
Conecte al servidor pomeo usando la siguiente dirección IP: 200.54.220.141 puerto 22.
Para Windows puede conectarse con el software PuTTy y las claves de acceso dadas en la primera clase.
2 - Configurar bioconda e instalar programas para análisis
Para configurar el canal bioconda se debe ejecutar el siguiente comando
conda config --add channels biocondaInstalar los software bwa
conda install -c bioconda bwa
Instalar samtools requiere un proceso de instalación un poco más largo, siga cada paso con atención.
conda install -c bioconda samtools
conda config --add channels bioconda
conda config --add channels conda-forge
conda install samtools==1.113 - Verificar directorios de instalación
Si desea constatar la ruta de instalación de un programa determinado puedes colocar el comando “whereis” en la terminal más el programa del que desea obtener la ruta.
whereis sratoolkit
whereis samtools
whereis bwa La salida de cada comando indicará la ruta de instalación como sigue, para cada programa del que desees obtener información.
/home2/usuario/miniconda3/bin/bwa
4 - Creación de directorio de trabajo y descarga de datos para alineamiento
Considerando que nuestras secuencias fastq originales tuvieron muy buena calidad, puedes trabajar directamente con ellas. Para ambos casos debes trasladar los archivos a una nueva carpeta, así:
Crear una carpeta denominada en tu usuario de home2
mkdir alineamientoIngresa a tu carpeta “alineamiento” y transfiere los archivos de clase anterior colocando los siguientes comandos en la terminal:
mv /home2/usuario/SRA_samples/SRR2006763/SRR2006763_1.fastq /home2/usuario/alineamiento/
mv /home2/usuario/SRA_samples/SRR2006763/SRR2006763_2.fastq /home2/usuario/alineamiento/Lista tu carpeta de alineamiento para verificar que tienes lo necesario para el alineamiento, hasta ahora deben constar tus dos secuencias “SRR2006763_1.fastq” y “SRR2006763_2.fastq”
ls Ahora descargaremos el genoma de referencia de la mitocondria de Salmo salar en la misma carpeta de alineamiento, para ello:
En tu navegador entra al siguiente link
Encontrarás una tabla con el listado del genoma de referencia donde se incluyen todos los cromosomas y el genoma de la mitocondria y a cuyo identificardor RefSeq darás clic, como se ve en la siguiente imagen
Listado genoma Salmo salar
Dar clic en la opción FASTA localizada bajo el título e identificador RefSeq de la referencia
Fasta option
Enviar la secuencia FASTA del genoma mitocondrial a un archivo
Genoma mitocondrial
Ve al lugar de descargas de tu computador y encontrarás el archivo denominado “sequence.fasta” al que debes renombrar como “mt.fasta” y subirlo a POMEO, para ello puedes trabajar con WINSCP
Si trabajas con WINSCP debes iniciar sesión, guardar tus datos y conectarte como se indica en la siguiente imagen
WINSCP
Al conectar debes ingresar nuevamente tu contraseña y aparecerá la interfaz de tu servidor con las carpetas que tienes creadas, aquí ingresarás en tu carpeta de alineamiento y arrastrarás el archivo descagado del genoma mitocondrial a la misma.
Alternativamente puedes cargar el genoma iniciando sesión online en el servidor POMEO a traves del puerto 8787 precargado con Rstudio server.
http://200.54.220.141:8787/auth-sign-inCuando termines lo anterior puedes ingresar a POMEO, listar tu carpeta de alineamiento con ls y verás tu archivo “mt.fasta” junto con tus secuencias fastq.
5 - Indexación del genoma de referencia
Una vez que incluiste en tu carpeta de alineamiento todos los archivos descritos en pasos anteriores podemos proceder a la etapa inicial del alineamiento, que corresponde a la indexación del genoma de referencia con bwa usando el siguiente comando:
bwa index mt.fastaLa salida del comando dará como resultado 5 archivos con extensiones “amb”,“ann”,“bwt”,“pac” y “sa”
5 - Alineamiento
Para el alineamiento tendremos las siguientes etapas
Alineamiento de las secuencias contra el genoma de referencia, cuya salida será un archivo con extensión “.sam”
Conversión del archivo SAM a BAM
Inspeccionar el archivo .sam de salida
Ordenar lecturas alineadas por posición
Indexación con Samtools
Exploración de datos con Samtools
Para ejecutar todas las etapas anteriores en ese orden se debe crear un script con nano denominado aln_mt.sh
nano aln_mt.shEn el script coloca las siguientes instrucciones
#!/bin/bash -l
# para alinear tus dos secuencias fastq al genoma mitocondrial
bwa mem mt.fasta SRR2006763_1.fastq SRR2006763_2.fastq > SRR2006763.sam
# Transformar tu archivo sam a bam
samtools view -Sb -q 30 SRR2006763.sam > SRR2006763.bam
# ordenar tu archivo binario bam
samtools sort SRR2006763.bam -o SRR2006763.sort.bam
# indexar tu archivo bam
samtools index SRR2006763.sort.bam Ejecuta tu script con bash
bash aln_mt.shAhora que ya tienes tus archivos SAM/BAM puedes observar tu archivo sam con el comando less de linux (recuerda que es un archivo de texto plano)
less SRR2006763.samtambién puedes realizar un análisis estadístico estandar con los siguientes comandos
samtools flagstat SRR2006763.bam > muestra_stat.txtInvestiga otros comandos de samtools para explorar el alineamiento.
Algunas observaciones adicionales
BWA-MEM se recomienda para secuencias de más de ~70 pb, ya que, generalmente BWA-MEM es más tolerante con los errores dadas las secuencias más largas.
Si quieres ver lo que sucede en cada etapa de alineamiento, puedes colocar en la terminal cada comando del script de forma individual, e.g.
bwa mem mt.fasta SRR2006763_1.fastq SRR2006763_2.fastq > SRR2006763.sam# para alinear tus dos secuencias fastq al genoma mitocondrialsamtools view -Sb -q 30 SRR2006763.sam > SRR2006763.bam# Transformar tu archivo sam a bamsamtools sort SRR2006763.bam -o SRR2006763.sort.bam# ordenar tu archivo binario bamsamtools index SRR2006763.sort.bam# indexar tu archivo bamEl archivo bam no se puede visualizar directamente porque está comprimido y en binario.
Procura listar tu carpeta cada vez que realices un trabajo, una descarga, o ejecutes un script para comprobar tus salidas.
Si quieres comprobar los tamaños de tus archivos al mismo tiempo que los listas, puedes hacerlo con
ls -l -h
REFERENCIAS
Pham-Quoc, C., Kieu-Do, B., & Thinh, T. N. (2019). A high-performance FPGA-based BWA-MEM DNA sequence alignment. Wiley, 1-12.
Li, H., & Durbin, R. (2009). Fast and accurate short read alignment with Burrows-Wheeler Transform. Bioinformatics, 1754-60.
JHU. (2016). Johns Hopkins University. Obtenido de Manual of Bowtie 2 - Fast and sensitive read alignment: http://bowtie-bio.sourceforge.net/bowtie2/manual.shtml
